/*
 * Kotlin/Native N-Body Simulation
 * Fully corrected version
 */

import kotlin.math.sqrt
import kotlin.math.abs

class Body(
    var x: Double,
    var y: Double,
    var z: Double,
    var vx: Double,
    var vy: Double,
    var vz: Double,
    var mass: Double
)

object NBody {
    private const val PI = 3.141592653589793
    private val SOLAR_MASS = 4 * PI * PI
    private const val DAYS_PER_YEAR = 365.24

    private fun offsetMomentum(bodies: Array<Body>) {
        var px = 0.0
        var py = 0.0
        var pz = 0.0

        for (body in bodies) {
            px += body.vx * body.mass
            py += body.vy * body.mass
            pz += body.vz * body.mass
        }

        bodies[0].vx = -px / SOLAR_MASS
        bodies[0].vy = -py / SOLAR_MASS
        bodies[0].vz = -pz / SOLAR_MASS
    }

    private fun energy(bodies: Array<Body>): Double {
        var e = 0.0
        for (i in bodies.indices) {
            val bi = bodies[i]
            e += 0.5 * bi.mass * (bi.vx * bi.vx + bi.vy * bi.vy + bi.vz * bi.vz)

            for (j in i+1 until bodies.size) {
                val bj = bodies[j]
                val dx = bi.x - bj.x
                val dy = bi.y - bj.y
                val dz = bi.z - bj.z
                val distance = sqrt(dx * dx + dy * dy + dz * dz)
                e -= (bi.mass * bj.mass) / distance
            }
        }
        return e
    }

    private fun advance(bodies: Array<Body>, dt: Double) {
        for (i in bodies.indices) {
            for (j in i+1 until bodies.size) {
                val bi = bodies[i]
                val bj = bodies[j]
                val dx = bi.x - bj.x
                val dy = bi.y - bj.y
                val dz = bi.z - bj.z
                val dSquared = dx * dx + dy * dy + dz * dz
                val distance = sqrt(dSquared)
                val mag = dt / (dSquared * distance)

                bi.vx -= dx * bj.mass * mag
                bi.vy -= dy * bj.mass * mag
                bi.vz -= dz * bj.mass * mag

                bj.vx += dx * bi.mass * mag
                bj.vy += dy * bi.mass * mag
                bj.vz += dz * bi.mass * mag
            }
        }

        for (body in bodies) {
            body.x += dt * body.vx
            body.y += dt * body.vy
            body.z += dt * body.vz
        }
    }

    private fun printFormatted(value: Double) {
        val integerPart = value.toLong()
        val fractional = abs((value - integerPart) * 1e9).toLong()
        val fractionalStr = fractional.toString().padStart(9, '0')
        val sign = if (value < 0) "-" else ""
        println("$sign${abs(integerPart)}.$fractionalStr")
    }

    fun run(args: Array<String>) {
        val bodies = arrayOf(
            Body(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, SOLAR_MASS), // Sun
            Body(
                4.84143144246472090e+00,
                -1.16032004402742839e+00,
                -1.03622044471123109e-01,
                1.66007664274403694e-03 * DAYS_PER_YEAR,
                7.69901118419740425e-03 * DAYS_PER_YEAR,
                -6.90460016972063023e-05 * DAYS_PER_YEAR,
                9.54791938424326609e-04 * SOLAR_MASS
            ), // Jupiter
            Body(
                8.34336671824457987e+00,
                4.12479856412430479e+00,
                -4.03523417114321381e-01,
                -2.76742510726862411e-03 * DAYS_PER_YEAR,
                4.99852801234917238e-03 * DAYS_PER_YEAR,
                2.30417297573763929e-05 * DAYS_PER_YEAR,
                2.85885980666130812e-04 * SOLAR_MASS
            ), // Saturn
            Body(
                1.28943695621391310e+01,
                -1.51111514016986312e+01,
                -2.23307578892655734e-01,
                2.96460137564761618e-03 * DAYS_PER_YEAR,
                2.37847173959480950e-03 * DAYS_PER_YEAR,
                -2.96589568540237556e-05 * DAYS_PER_YEAR,
                4.36624404335156298e-05 * SOLAR_MASS
            ), // Uranus
            Body(
                1.53796971148509165e+01,
                -2.59193146099879641e+01,
                1.79258772950371181e-01,
                2.68067772490389322e-03 * DAYS_PER_YEAR,
                1.62824170038242295e-03 * DAYS_PER_YEAR,
                -9.51592254519715870e-05 * DAYS_PER_YEAR,
                5.15138902046611451e-05 * SOLAR_MASS
            ) // Neptune
        )

        val n = args.getOrNull(0)?.toIntOrNull() ?: 1000

        offsetMomentum(bodies)
        printFormatted(energy(bodies))

        repeat(n) {
            advance(bodies, 0.01)
        }

        printFormatted(energy(bodies))
    }
}

fun main(args: Array<String>) {
    NBody.run(args)
}